#rm(list=ls())
#setwd("/Main Text-Figures/Figure5/")

#library(wnominate)
#library(pscl)

########################################

alpha <- .5
k <- .35
N <- 100
v <- 500
nsim <- 25

#eigen.all <- NULL
#eigen.rcv <- NULL
#sim.info.all <- NULL

for(j in 1:nsim){
  phi <- runif(1,.2,.8)
  
  b1 <- runif(v,0,1)
  b2 <- runif(v,0,1)
  sq1 <- runif(v,0,1)
  sq2 <- runif(v,0,1)
  
  npL <- round((1-phi)*N,digits=0)
  npR <- round((phi)*N,digits=0)
  
  d1 <- runif(1,.5,2)
  d2 <- runif(1,.5,2)
  
########################################
    
    hypo.pro.1d <- cbind(rep(1,length(b1)),rep(d1,length(b1)),b1,sq1)
    colnames(hypo.pro.1d) <- c("dim","d","b","sq")
    hypo.pro.2d <- cbind(rep(2,length(b2)),rep(d2,length(b2)),b2,sq2)
    colnames(hypo.pro.2d) <- c("dim","d","b","sq")
    hypo.pro.all <- as.data.frame(rbind(hypo.pro.1d,hypo.pro.2d))
    hypo.pro.all$id <- seq(1,nrow(hypo.pro.all),1)
    
########################################
    
    xL1 <- c(0,runif(npL-1,min=0-d1,max=d1))
    xL2 <- c(0,runif(npL-1,min=0-d2,max=d2))
    xR1 <- c(runif(npR-1,min=1-d1,max=1+d1),1) 
    xR2 <- c(runif(npR-1,min=1-d2,max=1+d2),1)
    
    left.party <- as.data.frame(cbind(xL1,xL2,rep(0,npL)))
    colnames(left.party) <- c("x1","x2","partyid")
    right.party <- as.data.frame(cbind(xR1,xR2,rep(1,npR)))
    colnames(right.party) <- c("x1","x2","partyid")
    all.leg <- as.data.frame(rbind(left.party,right.party))
    
########################################
#votes?

    vote.func1 <- function(x1,b1,sq1){
      as.numeric(abs(x1-b1)  < abs(x1-sq1) )
    }
    
    vote.func2 <- function(x2,b2,sq2){
      as.numeric(abs(x2-b2)  < abs(x2-sq2) )
    }
    
    out <- matrix(data=NA,nrow=N,ncol=v*2);dim(out)
    for(i in 1:N){
      for(j in 1:v){
        out[i,j] <- vote.func1(all.leg$x1[i],b1[j],sq1[j])
      }
    }
    for(i in 1:N){
      for(j in (v+1):(2*v)){
        out[i,j] <- vote.func2(all.leg$x2[i],b2[j-v],sq2[j-v])
      }
    }
########################################
#pass?
    
    passes <- as.numeric(colSums(out) > .5*N)
    hypo.pro.all <- as.data.frame(cbind(hypo.pro.all,passes))
    for (i in 1:nrow(hypo.pro.all)){
      hypo.pro.all$outcome[i] <- if(hypo.pro.all$passes[i]==0) hypo.pro.all$sq[i] else hypo.pro.all$b[i]
    }
    
########################################

    left.aligned <- rep(NA,length=ncol(out))
    for(i in 1:ncol(out)){
      left.aligned[i] <- round((sum(as.numeric(out[2:npL,i]==out[1,i])))/npL,4)
    }
    right.aligned <- rep(NA,length=ncol(out))
    for(i in 1:ncol(out)){
      right.aligned[i] <- round((sum(as.numeric(out[(npL+1):(N-1),i]==out[N,i])))/npR,4)
    }
    
    out.cohesion <- matrix(data=NA,nrow=N,ncol=v*2)
    for(j in 1:(v*2)){
      for(i in 1:npL){
        out.cohesion[i,j] <- if (out[i,j] == out[1,j]) {
          left.aligned[j] - right.aligned[j]
        } else {
          0
        }
      }
    }
    for(j in 1:(v*2)){
      for(i in (npL+1):N){
        out.cohesion[i,j] <- if (out[i,j] == out[N,j]) {
          right.aligned[j]-left.aligned[j]
        } else {
          0
        }
      }
    }
    
    out.cohesion.abs <- matrix(data=NA,nrow=N,ncol=v*2)
    for(j in 1:(v*2)){
      for(i in 1:npL){
        out.cohesion.abs[i,j] <- if (out[i,j] == out[1,j]) {
          left.aligned[j] 
        } else {
          0
        }
      }
    }
    for(j in 1:(v*2)){
      for(i in (npL+1):N){
        out.cohesion.abs[i,j] <- if (out[i,j] == out[N,j]) {
          right.aligned[j]
        } else {
          0
        }
      }
    }

########################################
#RCV requested?
    
    out.cohesion.kmat <- matrix(data=NA,nrow=N,ncol=v*2)
    for(i in 1:N){
      for(j in 1:(v*2)){
        out.cohesion.kmat[i,j] <- as.numeric(out.cohesion[i,j] > k)
      }
    }
    hypo.pro.all$left.cscore <- out.cohesion[1,] 
    hypo.pro.all$right.cscore <- out.cohesion[N,]
    
    hypo.pro.all$left.cscore.abs <- out.cohesion.abs[1,] 
    hypo.pro.all$right.cscore.abs <- out.cohesion.abs[N,]
    
    hypo.pro.all$left.rcv.request <- out.cohesion.kmat[1,]
    hypo.pro.all$right.rcv.request <- out.cohesion.kmat[N,]
    
    hypo.pro.all$rcv.request <- as.numeric(colSums(out.cohesion.kmat)>0)

########################################
#proposed?
    
    for(i in 1:nrow(hypo.pro.all)){
      hypo.pro.all$proposing.party[i] <- if(hypo.pro.all$b[i] < hypo.pro.all$sq[i]) 0 else 1 
    }
    
    hypo.pro.all$proposed.left <- rep(0,length(nrow(hypo.pro.all)))
    for(i in 1:nrow(hypo.pro.all)){
      if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$left.rcv.request[i]==1 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$b[i])+hypo.pro.all$left.cscore[i]-k > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$left.rcv.request[i]==1 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$sq[i])+hypo.pro.all$left.cscore[i]-k > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$rcv.request[i]==0 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$b[i]) > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$rcv.request[i]==0 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$sq[i]) > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$right.rcv.request[i]==1 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$b[i])+hypo.pro.all$left.cscore[i] > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==0 & hypo.pro.all$right.rcv.request[i]==1 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.left[i] <- as.numeric(alpha*-abs(0-hypo.pro.all$sq[i])+hypo.pro.all$left.cscore[i] > alpha*-abs(0-hypo.pro.all$sq[i]))
      } else {
        0
      }
    }
    
    hypo.pro.all$proposed.right <- rep(0,length(nrow(hypo.pro.all)))
    for(i in 1:nrow(hypo.pro.all)){
      if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$right.rcv.request[i]==1 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$b[i])+hypo.pro.all$right.cscore[i]-k > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$right.rcv.request[i]==1 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$sq[i])+hypo.pro.all$right.cscore[i]-k > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$rcv.request[i]==0 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$b[i]) > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$rcv.request[i]==0 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$sq[i]) > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$left.rcv.request[i]==1 & hypo.pro.all$passes[i]==1) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$b[i])+hypo.pro.all$right.cscore[i] > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else if (hypo.pro.all$proposing.party[i]==1 & hypo.pro.all$left.rcv.request[i]==1 & hypo.pro.all$passes[i]==0) {
        hypo.pro.all$proposed.right[i] <- as.numeric(alpha*-abs(1-hypo.pro.all$sq[i])+hypo.pro.all$right.cscore[i] > alpha*-abs(1-hypo.pro.all$sq[i]))
      } else {
        0
      }
    }
    
    hypo.pro.all$proposed <- hypo.pro.all$proposed.left+hypo.pro.all$proposed.right
    
    all.proposals <- hypo.pro.all[which(hypo.pro.all$proposed==1),]
      all.proposals.left <- all.proposals[which(all.proposals$proposed.left==1),]
      all.proposals.right <- all.proposals[which(all.proposals$proposed.right==1),]
    rcv.proposals <- hypo.pro.all[which(hypo.pro.all$proposed==1 & hypo.pro.all$rcv.request==1),]
      rcv.proposals.left <- rcv.proposals[which(rcv.proposals$proposed.left==1),]
      rcv.proposals.right <- rcv.proposals[which(rcv.proposals$proposed.right==1),]

    vote.matrix.all <- out[,all.proposals$id]
      vote.matrix.all.left <- out[,all.proposals.left$id]
      vote.matrix.all.right <- out[,all.proposals.right$id]
    vote.matrix.rcv <- out[,rcv.proposals$id]
      vote.matrix.rcv.left <- out[,rcv.proposals.left$id]
      vote.matrix.rcv.right <- out[,rcv.proposals.right$id]
    
########################################
#Dimensionality
    
    votemat.all <- rollcall(vote.matrix.all, yea=1, nay=0)
    votemat.rcv <- rollcall(vote.matrix.rcv, yea=1, nay=0)
    
    res.all <- wnominate(votemat.all,polarity = c(100,100),dims=2,
                         verbose=FALSE)
    res.rcv <- wnominate(votemat.rcv,polarity = c(100,100),dims=2,
                         verbose=FALSE)
  
    eigen.all <- cbind(eigen.all,res.all$eigenvalues)
    eigen.rcv <- cbind(eigen.rcv,res.rcv$eigenvalues)
    
    temp.all <- as.numeric(table(all.proposals$d))
    temp.rcv <- as.numeric(table(rcv.proposals$d))
    sim.info <- c(nrow(all.proposals),nrow(rcv.proposals),phi,d1,d2,temp.all[2],temp.all[1],temp.rcv[2],temp.rcv[1])
    sim.info.all <- rbind(sim.info.all,sim.info)

}

eigen.all <- eigen.all[,1:25]
eigen.rcv <- eigen.rcv[,1:25]

colnames(sim.info.all) <- c("non.rcv","rcv","phi","d1","d2","all.votes.d1","all.votes.d2","rcv.votes.d1","rcv.votes.d2")
sim.info.all <- as.data.frame(sim.info.all)
sim.info.all <- sim.info.all[1:25,]
sim.info.all$total <- sim.info.all$non.rcv+sim.info.all$rcv
sim.info.all$percent.rcv <- (sim.info.all$rcv / sim.info.all$total)*100

load("dimensionality-fig5-data.RData")

pdf("Fig5_.pdf",width=8,height=4)
par(mfrow=c(1,2),mar=c(4,4,4,1))
plot(eigen.all[,1],xlim=c(1,10),ylim=c(0,20),type='l',ylab="",xlab="",yaxt='n',col="grey60",lwd=1,xaxt='n')
abline(v=2,col="black",lwd=4)
for(i in 1:ncol(eigen.all)){
  lines(eigen.all[,i],col="grey60",lwd=1)
}
axis(side=2,at=seq(0,20,5),labels=seq(0,20,5),las=1)
mtext(side=3,line=2.25,"(a) All Votes",cex=1.5)
mtext(side=3,line=.75,print(paste("Average Total Votes  = ",round(mean(sim.info.all$total),0))),cex=1)
text(x=10,y=19,print(paste("k = ",k)),cex=1,pos=2)
text(x=10,y=15,print(paste("Missed Second\n Dimension: ","0-4%" )),pos=2)
mtext(side=2,line=3,"Eigenvalues")
mtext(side=1,line=2.5,"Number of Dimensions")
axis(side=1,at=seq(2,10,2),labels=seq(2,10,2))
axis(side=1,at=c(1,3),labels=c(1,3))
abline(h=1,lty=2,col="black")

par(mar=c(4,1,4,4))
plot(eigen.rcv[,1],xlim=c(1,10),pch=16,ylim=c(0,20),type='l',ylab="",xlab="",yaxt='n',col="grey60",lwd=1,xaxt='n')
abline(v=2,col="black",lwd=4)
for(i in 1:ncol(eigen.rcv)){
  lines(eigen.rcv[,i],col="grey60",lwd=1)
}
axis(side=4,at=seq(0,20,5),labels=seq(0,20,5),las=1)
mtext(side=3,line=2.25,"(b) Observed Roll Call Votes",cex=1.5)
mtext(side=3,line=.75,print(paste("Average % RCV  = ",round(mean(sim.info.all$percent.rcv),0),"%")),cex=1)
text(x=10,y=19,print(paste("k = ",k)),cex=1,pos=2)
text(x=10,y=15,print(paste0("Missed Second\n Dimension: ",((as.numeric(table(eigen.rcv[2,]<=1)[2])/25)*100),"-44%" )),pos=2)
mtext(side=4,line=3,"Eigenvalues")
mtext(side=1,line=2.5,"Number of Dimensions")
axis(side=1,at=seq(2,10,2),labels=seq(2,10,2))
axis(side=1,at=c(1,3),labels=c(1,3))
abline(h=1,lty=2,col="black")
dev.off()


########################################

mean(sim.info.all$rcv)

mean(sim.info.all$rcv / (sim.info.all$non.rcv+sim.info.all$rcv))
(sim.info.all$rcv / (sim.info.all$non.rcv+sim.info.all$rcv))

table(eigen.rcv[2,]<1)
table(eigen.all[2,]<1)

table(eigen.rcv < 2)
table(eigen.all < 2)

summary(sim.info.all$all.votes.d1)
summary(sim.info.all$all.votes.d1 / (sim.info.all$all.votes.d1 + sim.info.all$all.votes.d2))
summary(sim.info.all$rcv.votes.d1)
summary(sim.info.all$rcv.votes.d2)

########################################
#votes coming from what dimension?

sim.info.all2 <- sim.info.all

for(i in 1:nrow(sim.info.all2)){
sim.info.all2$all.votes.d1b[i] <- if (sim.info.all2$d1[i] > sim.info.all2$d2[i]) {
  sim.info.all2$all.votes.d1[i]
} else {
  sim.info.all2$all.votes.d2[i]
}
}

for(i in 1:nrow(sim.info.all2)){
  sim.info.all2$all.votes.d2b[i] <- if (sim.info.all2$d1[i] > sim.info.all2$d2[i]) {
  sim.info.all2$all.votes.d2[i]
} else {
  sim.info.all2$all.votes.d1[i]
}
}

for(i in 1:nrow(sim.info.all2)){
  sim.info.all2$rcv.votes.d1b[i] <- if (sim.info.all2$d1[i] > sim.info.all2$d2[i]) {
    sim.info.all2$rcv.votes.d1[i]
  } else {
    sim.info.all2$rcv.votes.d2[i]
  }
}
for(i in 1:nrow(sim.info.all2)){
  sim.info.all2$rcv.votes.d2b[i] <- if (sim.info.all2$d1[i] > sim.info.all2$d2[i]) {
    sim.info.all2$rcv.votes.d2[i]
  } else {
    sim.info.all2$rcv.votes.d1[i]
  }
}

sim.info.all2$percent.rcv.d1 <- sim.info.all2$rcv.votes.d1/(sim.info.all2$rcv.votes.d1 + sim.info.all2$rcv.votes.d2)
mean(sim.info.all2$percent.rcv.d1,na.rm=TRUE)
sim.info.all2$percent.rcv.d2 <- sim.info.all2$rcv.votes.d2/(sim.info.all2$rcv.votes.d1 + sim.info.all2$rcv.votes.d2)
mean(sim.info.all2$percent.rcv.d2,na.rm=TRUE)

sim.info.all2$percent.all.d1 <- sim.info.all2$all.votes.d1/(sim.info.all2$all.votes.d1 + sim.info.all2$all.votes.d2)
mean(sim.info.all2$percent.all.d1,na.rm=TRUE)
sim.info.all2$percent.all.d2 <- sim.info.all2$all.votes.d2/(sim.info.all2$all.votes.d1 + sim.info.all2$all.votes.d2)
mean(sim.info.all2$percent.all.d2,na.rm=TRUE)

mean(sim.info.all$total)
mean(sim.info.all$rcv)
mean(sim.info.all$rcv)/mean(sim.info.all$total)

